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A practical review of the Adomian 
decomposition method: computer 
implementation aspects 


A. Molabahrami 


Abstract 


In this paper, a practical review of the Adomian decomposition method, 
to extend the procedure to handle the strongly nonlinear problems under the 
mixed conditions, is given and the convergence of the algorithm is proved. 
For this respect, a new and simple way to generate the Adomian polynomials, 
for a general nonlinear function, is proposed. The proposed procedure, pro- 
vides an explicit formula to calculate the Adomian polynomials of a nonlinear 
function. The efficiency of the approach will be shown by applying the pro- 
cedure on several interesting integro-differential problems. The Mathematica 
programs generating the Adomian polynomials and Adomian solutions based 
on the procedures in this paper are designed. 


Keywords: Adomian decomposition method; Adomian polynomials; Non- 
linear integro-differential problems; Series solution; Strongly nonlinear prob- 
lems; Explicit machine computation and programs. 


1 Introduction 


To construct series pattern solution for a problem with strong nonlinearity, 
it is necessary to construct nonlinear terms of the governing equation in 
the form of a series by using the components of the solution series. To the 
end, one of the best and suitable way is to use the Adomian polynomials. 
The so-called Adomian polynomials are used to deduce the recursive relation 
during the implementation of the Adomian decomposition method (ADM) 
while solving nonlinear problems. The main aim of the present paper is to 
provide a simple and new method to handle a strongly nonlinear problem 
by using ADM in the frame of a symbolic computer program so that by 
giving linear operator, it generates: initial guess, integral inverse of the linear 
operator, recursive relation and the terms of solution series automatically. 
To achieve this purpose, we first propose an explicit formula to calculate the 
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Adomian polynomials and the Adomian series of a general nonlinear function 
and implement the proposed algorithms in Mathematica. In this respect, we 
first outline the modifications of some definitions as already given in [8]. Let 
u be a function of the parameter A, whose Maclaurin series is given by 


+00 
u(A) = Un”. (1) 
n=0 


This series is called the parametric series of u. Let @ be a function of the 
parameter A, the mth-order parametric derivative of ¢ is 


1 O"¢ 


Dml[¢] = mi ave 


(2) 


A=0 


where m > 0 is an integer. The mth-order Adomian polynomial of ¢ is 


Am(9(u)) = Dm [6(u(A))], (3) 


where m > 0 is an integer and A,,(¢(u)) = Am(b(u); Ug, U1, ---) Um). 


Remark 1. For the case 0 <  < 1, the parametric series (1) and parametric 
derivative (2) reduce to homotopy series and homotopy derivative respectively 
[8]. 


Several algorithms [3,5,12,13] for symbolic programming have since been 
devised to efficiently generate the Adomian polynomials quickly to high or- 
ders, for example, a convenient formula for the Adomian polynomials is the 
rule of Rach, which reads (see Page 16 in [1] and Page 51 in [2]) 


n 


k=1 


where the C(k,n) are the sums of all possible products of & components 
uj,t = 1,2,.....-—k+1, whose subscripts sum to n, divided by the factorial 
of the number of repeated subscripts. An equivalent expression of Equation 


An(f(u)) = S- pf tee P Gis) II a n>1. (5) 
s=1 *S° 


pit2p2t+...+npn=n 


In the present paper, we propose an explicit formula to calculate the 
C(k,n) in (4) and we show that for rapid computer-generation of the Ado- 
mian polynomials there is no need to use the (5). 
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2 Adomian polynomials for a general function 

In this section, we first outline two theorems as already given in [7,10]. Then 
using them, we propose a new theorem which provides a new and simple 


way to calculate the Adomian polynomials for a general smooth function. 
Here, we mention them with a minor modification. For simplicity, we use the 


following notation 
m 
Umn = S- ur. 
i=n 


Theorem 1. For function f(u) =u", the corresponding mth-order Adomian 
polynomial is given by 


m TL r2 
Am (u ) — Um—r1 Ur i—re Urg—rg +e 


r1=0 r2=0 r3=0 
Tk—-3 Tk-2 
5 Urp—3—Tr-2 5 Urp—o—re—-1Urp_13 (6) 
TR-2=0 TR-1=0 


where m > 0 and k > 0 are positive integers. 


Proof. Considering the definition (2), refer to [10]. 


Theorem 2. For parametric series (1), it holds 


Dm [f(u(A))] = Dm [F(Um,o)] 5 


where f is a smooth function. 


Proof. Refer to [7]. 


Corollary 1. From Theorem 1, we find 


+00 h +00 
u®(Xr) = (s wa") =u + $7 Am(u*)A”, (7) 
n=0 m=1 
where the Adomian polynomials Am(u") are given by (6) . 


Remark 2. It is clear that A,,(u*) in Theorem 1 can easily be calculated 
by a simple code by using a symbolic software such as Mathematica. For this 
respect, the Code. 1, reported in Appendix, can be used in Mathematica. 
For instance, by ADPforPowerLaw[4,6], the A4(u®) is calculated as follows 


Ag(u®) = 15ugut + 60uguzus + 15ugu3 + 30upurus + 6ugua. 


Corollary 2. From Theorem 2, form >n, we find 
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Dm [(f (too,n)] =Dmn [f (@m,n)] ) 


where f is a smooth function. 
Corollary 3. From Corollary 2 and Theorem 1, form > k, we find 


m-1l m1 Tk-3 Tk—-2 k-2 
s k 
Dm [(@m,2) = S S Pat S S Urn I] Urj;—rj419 (8) 
r1=0 r2=0 Te—-2=0 Tk—-1=0 j=0 
r2AT1 Te-2ATK—-3 TR-1FATk—-2 


where ro = ™. 
Corollary 4. It is easy to see that 


0, m<k, 


Drs [(Gma)*] = {Ome (9) 


Corollary 5. Let kn >m+1 andn>1, we find 


Din ey =, 


Remark 3. A sample Mathematica program for D,,  (im,1)*) in (8) is 


given by the Code. 2 reported in Appendix. An alternative way is to use 
the Theorem 1 by taking ug = 0. To achieve this purpose, In Code 1, in the 
last command, the Expand[Dz,m] is replaced by Expand|Dz,m]/.uo > 0. 

The following theorem provides a suitable and simple way to construct a 
recurrent relation for a general smooth function appeared within a structure 
with nonlinear terms in the equations. 


Theorem 3. Assume that f(u) has the Taylor expansion with respect to 
uo, then 


mek) (a) 
Am (F(a) =), [@ina)*) (10) 


ra k) 


Proof. Expanding f(u) in Taylor series with respect to uo, one has 


+20 #(k) (4 
F(t) = F(to) +S EO) (wy — uy. (11) 
k=1 
From the (11), we have 
Arn LF (w)] =D [> EM (wr) — ra) 
k=1 : 


recalling the Corollaries 5 and 2, we find 
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HON i aa (OVEN . 
Am (F (a) = ,, [(ua) ~ 10)"] = 2 MD, [Gna]: 
k=1 


k=1 


This ends the proof. 


Remark 4. The expansion (11) is rigorously guaranteed by the Cauchy- 
Kovaleski theorem. In fact, by using this theorem the Taylor expansion 
about the function uo(#) can be seen as an expansion in Banach space [4]. 


Remark 5. The expression D,,, (m1) in (10) can easily be calculated by 
(8) or Theorem 1 with taking up = 0. 


Corollary 5. From Theorem 3, we find 


+00 
F(u(A)) = f(uo) + D5 wd”, (12) 


where Wm = Am(f(u)) can be calculated by (10). 


Remark 6. To calculate Ay, (f(u)), the Code. 3, which is given in the Ap- 
pendix, can be used in Mathematica. Using the AdomianPolynomial|[f|u], m], 
form = 1, 2,3, 4, the corresponding Adomian polynomials are given as follows 


Ai(f(u)) = uf’ (uo), 
Ao(f(u)) = u2f’(uo) + sud f” (uo), 


A3(f(u)) = us f’(uo) + uiu2f” (uo) + suff!” (uo), 


Aa(f(u)) = ua f’ (uo) + #4 (24u3 + 48uiu3) f’”(uo) + dutusf’” (uo) + quik (uo). 


Also, the Code. 4, reported in Appendix, shows an alternative way to calcu- 
late the Adomian polynomials by using the (3) and Corollary 2. 


2.1 An improvement for the Rach’s rule 


According to (4) and Theorem 3, we find 


C(k,n) = <D, [(na)*) 


thus, the Rach’s rule gets its explicit presentation. 
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3 A practical review of the Adomian decomposition 
method 


In this section, we propose a practical review of the ADM and implementation 
it as an automatic program in the frame of a symbolic software such as 
Mathematica. 

Adomian decomposition method [1], was first proposed by Adomian in 
1988 and was further developed and improved by Adomian [2,3]. To illustrate 
the ADM for solving a general nonlinear problem, consider the following 


nonlinear problem 
N(u) =f, 
(Bote _ 


where N is a general nonlinear operator, B is a linear initial/boundary op- 
erator, u is the unknown function that will be determined and f and g are 
given functions. An especial case of (13), f = 0, was given in [9]. The ADM 
consists in looking for the solution of Equation (13) in the series form 


+00 
w= t+ Stn (14) 
n=1 


where up is an initial guess and has the property 


(uo) = fo, (15) 


and the other components of the solution series (14) have the property 


B(un) =0,n> 1. (16) 


Thus, from (15) and (16) the solution series given by (14), satisfies the condi- 
tions of (13). For the nonlinear conditions, the ADM needs an improvement. 

To construct series pattern solution, (14) by ADM, in the first step the 
nonlinear operator NV is decomposed as 


N(u) = L(u) + N(u), (17) 


where L is a linear operator and N = N — L. To continue in ADM, the 
nonlinear operator N is decomposed as 


+00 
N(u) => An, (18) 
n=0 


where A,, is the nth-order Adomian polynomial of N. Using (17), the Equa- 
tion (13) becomes 


L(u) + N(u) = f, 
{ (u) = fo, 9) 
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by ignoring the condition B(u) = fo, the solution of Equation (19) satisfies 


u= Ug + L"[f] - L"[N(u)], (20) 


where w, is the general solution of the linear equation L(u) = 0 and L~ is 
the integral inverse of L. Now, by substituting (14) in (20) and considering 
(18), we get 


Uo = Ug; 
le = (1=xn)L7'[f]- £7 An-i], 22 1. (21) 
where 
0, n<1, 
a 1, n>2. 


Recalling the condition B(u) = fo, we have 


ee B(ug) = ux, (22) 


B(un) =0, n>1. 


Thus, the recurrent relation to find a series pattern solution of Equation (13) 
by means of ADM is 


uo =U", 
he gait, 1], Btn) = 0, 23) 


let uParticular = L~*[—Ap_1] + (1 — xn) Z7[f], then, for n > 1, we find 


Up = ur", 
(24) 
Un = UuParticular + uParticular Star, 


where uParticular Star = B (L~*[An-1])—(1 — Xn) (B(L™1[f])). It is impor- 
tant to notice that the definition of the integral inverse L~! in (23), depends 
on the condition (16). L~' in (24) can be defined such that the L~"[x] gives 
a particular solution of the equation L(u) = *. 

Remark 7. It should be emphasized that the main step of the ADM is 


to choose a proper linear operator. According to (24), the linear operator 
should be defined so that 


1. The following problem has a solution 


L(uo) = 0, 
{ Ca (25) 


According to the Equations (20) and (21), the problem (25) gives the 
initial guess. 
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2. The following equation has a solution 
B(L~*[An—1]) — (1 — x(n)) B(L7"[f]) = Blu), (26) 


where n > 1. This equation helps to define the integral inverse of L. 


3. The solution series, given by (14), is convergent. For this respect, some 
conditions have been given in [11]. 


3.1 Decomposition series and convergence analysis 


An important hypothesis in ADM is the composition of the nonlinear operator 
N as shown in (18). There are some serious questions about (18) such as, 
from where does the series poke An derive? Does it always converge? Is 
its sum really N? What are the other series that we could use instead of 
ee A,? In order to answer these questions, and also, in order to explain 
some other issues, Gabet proposed a theory which explains and justifies the 
practical method [6]. In this subsection, we answer to the above mentioned 
questions by using the results of the section 2. 


To derive (18), from (12), under the assumptions of the Theorem 3, we 
have 


+00 
N(u(A)) =N(uo) + 55 And”, (27) 


where A,, is the nth-order Adomian polynomial of N. Let the parametric 


series (1) is convergent at \ = 1 and s = pa Un, then 
+00 
N(s) = Ap + 5) An. (28) 
n=1 


According to the (23), it is clear that if the solution series (14) is convergent 
to s, then the decomposition series ae A, converges to f — L(s). On 
the other hand, according to the (23), the decomposition series bra Ay, 
converges to N(s). Therefore, if the solution series given by ADM, of the 
form (14) generated by (24), is convergent, then it must be an exact solution 
of the considered nonlinear problem, denoted by (13). Bearing in mind the 
above discussion, we can express the following theorem on the convergency 
of the ADM. 


Theorem 4. Assume that the ADM solution series given by (14) and 
(24) is convergent, then, under the assumptions of the Theorem 8, it must be 
an exact solution of the (18). 
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4 Test examples 


To show the efficiency of Theorem 4, described in the previous section, 
some examples are presented. We consider the mth-order nonlinear integro- 
differential equation with variable coefficients 


m B 
So fr(a)u” (a) +2 i) K(a,t)F(u(t))dt = g(x), x € [a,b], (29) 


where F(u(a)) is a function of u(x), K(x,t) is the kernel of the integro- 
differential equation, is a parameter, g(x) is the data function, f,(x) is 
a function with respect to x, u(x) is the unknown function that will be 
determined. For the Fredholm and Volterra kinds we take G = b and 6 = x 
respectively. We consider Equation (29) under the mixed conditions 


m-1 


S- (a,ju! (a) + bp jul” (b) + enjul”(c)) = pj, 7 =0,1,2,..,m—1, (30) 
r=0 


where a,j, bj, Cr; and 1; are constants, and c, a < c < b, is a constant. For 
the linear case, it is assumed that F(u(#)) = u(x). The computations will 
be performed using the program ADMforIDEs reported in the appendix. 
The program ADMforIDEs has designed in a general manner so that, for 
a given problem and its related linear operator, it calculates the integral 
inverse of the linear operator, initial guess, recursive relation and the terms 
of solution series automatically. 


For Examples 1-4, we choose the linear operator as follows 


Lu(2)) = “5u(a), (31) 


therefore, from the (22) and (31), we find 
uo(x) = a. (32) 


We emphasize that no attempt has been made here to obtain all solutions of 
a given problem. 


Example 1. Consider the following nonlinear problem 
(33) 


An exact solution is u(x) = x. Starting by (32), the recurrent relation (24) 
gives 
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Un(x) =0, n> 1, 


therefore, using the ADM, the exact solution of Equation (33) is obtained as 
follows 


u(a) = ug(a) + Ss Un(x) = @. (34) 


Computations have been carried out by program ADMforIDEs. The list 
of commands is as follows 

conditions|u_] := {(u[z]/.2 > 3) + (D[u[z], 2]/.2 + 0) - 3, 
2* (u[z]/.2 > 5) + (D[ulz], 2] /.2 > 1) — 2}; 

Leoef ficients = {0,0, 1}; 

Problemcoef ficients = {—1, x, 1}; 

ADMforIDEs[Sin|[u], x, x * t,0, 5,1, 2,5] 

Example 2. Consider the following nonlinear problem 


u" (x) + cul (x) — u(x) + 2x fy te“ Odt =a —2xe-™, 
(35) 
u(s) + u/(0) = 8, 2u($) + w’(1) =2. 


The exact solution is u(a) = x. Using the recurrent relation (24) with initial 
guess (32), we find 
Un(xz) =0, n> 1, 


therefore, the ADM gives the exact solution of Equation (35) as follows 


u(x) = uo(x) + S- Un(x) = @. (36) 


Computations have been carried out with the help of the program ADM- 
forIDEs. The list of commands is as follows 
conditions|[u_] := {(u[az]/.2 > 3) + (D[u[z], 2]/.2 + 0) — 3, 
2 (ule]/0 > 4) + (Dlule],e]/.2 + 1) - 2}; 

Leooef ficients = {0,0, 1}; 

Problemcoef ficients = {—1, x, 1}; 

ADMforIDEs[Exp[—u?], 2 — « * Exp[—2x?],2* a * t,0,2, 1, 2,5] 


Example 3. Consider the following nonlinear problem 


u’ (x) + vu' (x) — u(x) + (2m + 2)ax 2 iV imalaa C3 el aeeerceon 
(37) 
u(}) + ul(0) = 8, 2u() + wa) =2, 


where m is a positive integer. An exact solution is u(#) = x. Recalling the 
recurrent relation (24) and initial guess (32), we obtain 
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Un(xz) =0, n> 1, 


therefore, the ADM obtains the exact solution of Equation (37) as follows 


u(x) = uo(a) + os Un(x) = @. (38) 


To achieve the computations by program ADMforIDEs, the list of com- 
mands is as follows 

conditions|{u_] := {(u[x]/.2 > 3) + (D[u[z], 2]/.2 + 0) — 3, 
2 (uls|/x + 3) + (D[ulz],2]/.0 - 1) — 2}; 

Leoef ficients = {0,0, 1}; 

Problemcoef ficients = {—1, x, 1}; 

$Assumptions = m > 0&&Element|m, Integers]; 
Refine|[ADMforIDEs[u2™, x, (2m + 2) « x * t,0, 1,1, 2, 5]] 


Example 4. Consider the following nonlinear problem 
u(x) + vu’ (a) — u(x) 4 Amt yy f° ¢ 2x /u(t)dt =f, 


u(5) +w(0) = 5, 2u(Z) + w'(1) = 2, 


(39) 


where m is a positive integer. An exact solution is u(x) = x. Putting to use 
the recurrent relation (24) and initial guess (32), we get 


Un(z) =0, n> 1, 


therefore, the ADM leads to the exact solution of Equation (39) as follows 


u(x) = uo(a) + 2 Un(ax) = @. (40) 


For performing the computations with the help of the program ADM- 
forIDEs, the list of commands is as follows 
conditions|{u_] := {(u[x]/.2 > 3) + (D[u[z], 2]/.2 + 0) — 3, 
2 (uls|/x + 3) + (D[ulz],2]/.0 - 1) — 2}; 

Leoef ficients = {0,0, 1}; 

Problemcoef ficients = {—1, x, 1}; 

$Assumptions = m > 0&&Element|m, Integers]; 
Refine[ADMforIDEs | ?%/u, x, #@++ « x « t,0, 1, 1,2, 5]] 


2m 


Remark 8. In Examples 1-4, the first-order iteration gives the exact so- 
lution. 


Example 5. Consider the following problem 
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e2 
(41) 
u(0) = 1, 
we choose the linear operator as follows 
d 
L(u(z)) = = u(x) + ule), (42) 


therefore, from the (22) and (42), we find uo(x) = e~*. Putting to use the 
recurrent relation (24) and obtained initial guess, we have 


1 
Un(x) = AG = Lyte — Lee}, n> 1, 


therefore, the ADM leads to the exact solution of Equation (39) as follows 
1 
u(a) = u(x) + S° un(2) = aa (-e*+e™th+1+e). (43) 


For performing the computations with the help of the program ADM- 
forIDEs, the list of commands is as follows 

conditions|u_] := {(u[x]/.c + 0) — 1}; 

Leoef ficients = {1,1}; 

Problemcoef ficients = {1,1}; 

ADM forIDEs [u, 5 (2s — 1) ,1,0,1,-1,1,5] 


Remark 9. The program ADMforIDEs, reported in the Appendix, mon- 
itors some terms of the solution series, the number of the terms is given by 
SolutionOrder. To generate the general term of the solution series (if pos- 
sible) and the closed form of the solution series (if possible), the user must 
be added the following commands to the end of the list of commands . It is 
important to emphasis that the following commands must be used in Math- 
ematica 7 or later. 

ml =Input["By evaluating the obtained terms solution series, please input the index of 
the term in which it is possible to generate the general term of the solution series?"]; 
SolutionTerms = Join[Table[u[m], {m, m1, SolutionOrder}]]; 

SolutionCoefficents = FindSequenceFunction[SolutionTerms, n]J; 

Print ["un(x)=", TraditionalForm|[SolutionCoefficents]] ; 

FinalSolution = Sum[SolutionCoefficents, {n, m1, Infinity},GenerateConditions + True]; 
Soll = Simplify [ mt) un] + FinalSolution| ; 


n=0 


Print ["u(a)= a un(x)=", TraditionalForm|[Sol1]| 
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5 Concluding remarks 


By calculating the parametric derivative for a function of a more general type, 
described in Theorem 3, the Adomian decomposition method becomes a pow- 
erful analytic approach for obtaining convergent series solutions of strongly 
nonlinear problems governing physical models in applied science and engi- 
neering. Using the parametric derivative, some lemmas and theorems pro- 
vided in ADM papers in the literature have been unified and modified here 
via some new theorems and corollaries. Also, by means of the (24) the Ado- 
mian decomposition method was extended to handel the nonlinear problems 
with the mixed conditions. 
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App endix (Mathematica programs) 


Code. 1(A sample Mathematica program for Am(u") given by (6)) 


ADPforPowerLaw(|m_,k_] := Module[{},D1,j_ := uj; Di_jo = ub; 
J order 
D25_ = Uj—rUr; Di_,order_ c= Dy Uorder—rDi-1,r3 
r=0 r=0 


Print (AY, (“uk”) =", Expand[TraditionalForm[D,, ]]]]; 


Code. 2(A sample Mathematica program for Dy, [(am1)*] given by (8)) 


j-1 
DnonUO[m_, k_] := Module[{}, D1,;_ := If[j > 0, u;, 0); Di_ jo := 0; Daj_ = YO uj—rur; 
r=1 
order—1 
Di_ order = ye Uorder—rDi-1,r3 Print|” DY,” Cae oat) =", 


r=2 
Expand[TraditionalForm[D, m]]]]; 


Code. 3(A sample Mathematica program for Am(f(u)) given by (10)) 


AdomianPolynomial[function—,m_—] := Module[{}, D1,j_,2_ := If[j > 0,25, 0); 
g-1 order—1 

Di_,0,2_ = 0; Do,j_ je = > Lj-r@&r; De onder ae = > Lorder—rDi-1,7,23 
r=1 r=2 


HDoraer_ ?= If[m == 0, function/.u > uo, 
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psa ((D[ function, {u, n}]/.u 4 uo) /n! * Dr orderjul); 


n=1 


Print[ A‘ (", function,” ) =", Expand[TraditionalForm[H D»]]]]; 


Code. 4(A sample Mathematica program for Am(f(u)) by using the (3) & 
Corollary 2) 


AdomianPolynomial1|function_,m_] := Module[{}, glu—] := function; 
A[m1_] := Factorial (Diglul/-u — Sum[un * p”, {n, 0, m1}], {p, m1}]/.p — 0); 


Print[A’’ ,”” (", function,” ) =", Expand[TraditionalForm[A[m]]]]]; 


ADMforIDEs(A sample Mathematica program of ADM for solving IDEs 
given by (29)) 


ADM forIDEs|functionF_, functiong-, kernel_,a_, B-, A-, ProblemOrder_, SolutionOrder_| 
:= Module[{}, F[u_] := functionF; g[x_] := functiong; k[x_, t_] := kernel; 


(« — — — —Definition of the linear operator *) 
L[f] := Module[{}, Expand[S7P70hemOrder (Looe f ficients|[i + 1] « D[f, {x, }])]]; 

(« — — — —Definition of the integral inverse of linear operator — — — — — — — *) 
Linverse|f_] := Module[{ Linu, u, solution}, Linv = DSolve[L[u[z]] == f, ula], x]; 
solution = Linv|[1, 1, 2]]/.C[-] — 0; Expand[solution]]; 

Linverse[p_Plus] := Map[Linverse, p]; Linverse|c_* f_] := c * Linverse[f]/; FreeQ|c, 2]; 
(« — — — —Definition of the Adomian polynomials *) 


HDN{horder_] := Module[{}, D1,;_ := If[j > 0, u[y], 0]; Di_o := 0; 
D2;_ := If[j > 0, vith ulj — 7] * ulr]]; Di_m_ = pias ulm — 7] * Di-1,r3 


DH _, function. ‘= Uf[horder == 0, function/.u — uO], 
Expand[}7""_, ((D[ function, {u, n}]/.u > u[0]) /(n!) * Drm) III; 


n=1 
* — — — —Definition of the initial guess *) 
ugStar = DSolve[L[u[x]] == 0, conditions[u] == 0, u[z], 2]; 


u[0] = First[u[a]/.ugStar]; 

x — — — —Definition of the ug x) 
ug = DSolve[L[u[x]] == 0, ula], x]; 

ugl = First[u[2]/.ug]; 


* — — — —Some needed commands *) 
X[m_] := If[m < 1,0, 1];C1 = Arrary[C, ProblemOrder]; w[x_] := ug]; 
Which|ProblemOrder == 1, constants[u_] := First[conditions[u]]/.[|— > 0, 
ProblemOrder > 1, constants[u_] := conditions[u]/.-[]— > 0]; 


fll = constants[ul] — conditions[w]; 

(« — — — —Main block *) 
For[m = 1,m < SolutionOrder,m ++, HDN{[m — 1; 

{A[m-1] = ye (Problemcoef ficients|[i + 1]] — Lcoef ficients|[i + 1]]) 

* D[ulm — 1], {x, }] + A * Integrate[k[a, t] * ((DHm_1,F[u])/-e > t), {t,4, Bh], 
uParticular = Linverse[—A[m — 1] + (1 — X[m]) * g[z]], wl[x_] := uParticular, 

f12 = conditions[w1] — constants[ul], ug2 = Solve[{ f11 == f12}, C1], 
uParticular Star = First[ug1//.ug2], 
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ulm] = uParticular + uParticularStar, 
Print!’ m ,/’ =", TraditionalForm[Collect[u[m], x]]]}]]; 
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